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Abstract 

A theoretical framework for the calculation of shear and bulk vis- 
cosities of hadronic matter at finite temperature is presented. The 
framework is based on the quasi-particle picture. It allows for an arbi- 
trary number of hadron species with point-like interactions, and allows 
for both elastic and inelastic collisions. Detailed balance is ensured. 
The particles have temperature dependent masses arising from mean 
field or potential effects, which maintains self-consistency between the 
equation of state and the transport coefficients. As an example, we 
calculate the shear and bulk viscosity in the linear a model. The ratio 
of shear viscosity to entropy density shows a minimum in the vicin- 
ity of a rapid crossover transition, while the ratio of bulk viscosity to 
entropy density shows a maximum. 

PACS numbers: ll.lO.Wx, 25.75.Nq, 51.20.+d, 25.75.-q 

Keywords: Shear viscosity, bulk viscosity, hadronic matter at finite temper- 
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1 Introduction 



One of the amazing experimental discoveries of measurements of heavy ion 
colUsions at the Relativistic Heavy Ion Colhder (RHIC) is the surprising 
amount of coUective flow exhibited by the outgoing hadrons. Collective flow 
is observed in both the single-particle transverse momentum distribution [1] 
(radial flow) and in the asymmetric azimuthal distribution around the beam 
axis [2] (elliptic flow). It is now generally accepted that collective flow is 
mostly generated early in the nucleus-nucleus collision and is present before 
partons fragment or coalesce into hadrons [3]. The quark-gluon matter cre- 
ated in these collisions must be strongly interacting, unlike the type of weakly 
interacting quark-gluon plasma expected to occur at very high temperatures 
on the basis of asymptotic freedom [4]. Perfect fluid dynamics with zero 
viscosity reproduces the measurements of radial and elliptic flow quite well 
up to transverse momenta on the order 1.5 GeV/c [5]. These results have 
been interpreted as strong indicators of early thermalization and collective 
flow on a time scale of several fm/c. 

An amazing theoretical discovery was made by Kovtun, Son and Starinets 
[6]. They showed that certain special field theories (AdS/CFT or Anti- 
deSitter/Conformal Field Theory) that are dual to black branes in higher 
space-time dimensions [7]- [9] have the ratio of shear viscosity to entropy 
density rj/s = I/Ati (in natural units with ft — — c — 1). The connec- 
tion between transport coefficients and gravity arises because both involve 
commutators of the stress-energy-momentum tensor. They conjectured that 
all substances have this value as a lower limit, and gave as examples various 
atomic and molecular systems. In fact, it had been argued much earlier that 
any substance should a lower bound on rj/s because of the uncertainty prin- 
ciple [10]. Is the RHIC data telling us that the created matter has a very 
small viscosity, the minimal value of r]/s, that it is a perfect fluid! 

The relatively good agreement between perfect fluid calculations and ex- 
perimental data for hadrons of low to medium transverse momentum at RHIC 
suggests that the viscosity is small. However, it cannot be zero. Indeed, cal- 
culations within AdS/CFT suggest that rj > s/in [11]-[17]. (Whether or 
not this is a rigorous lower bound is still an open question [18]-[20].) There 
are strong theoretical arguments, and evidence from atomic and molecular 
systems, that r]/s should be a minimum in the vicinity of the phase tran- 
sition or rapid crossover between hadronic matter and quark-gluon plasma 
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[21, 22], and that the ratio of bulk viscosity to entropy density (/s should be 
a maximum there [23]. See also [24]-[27]. 

It ought to be possible to extract numerical values of the viscosities in 
heavy ion collisions via scaling violations to perfect fluid flow predictions 
[28]- [31]. The program is to solve relativistic viscous fluid equations, with ap- 
propriate initial conditions and with a hadron cascade afterburner [32], over 
a range of beam energies and nuclei and extract ri{T)/s{T) and C,{T)/s{T) 
from comparison with data. Thus, sufficiently precise calculations and mea- 
surements should allow for a determination of the ratio r]/s as well as the 
ratio of bulk viscosity to entropy density C/s as functions of temperature, 
and that these ratios can pinpoint the location of the phase transition or 
rapid crossover from hadronic to quark and gluon matter. This is a different 
method than trying to infer the equation of state of QCD in the form of 
pressure P as a function of temperature T or energy density e. Because of 
advances in both theory and computation, vigorous activities are currently 
underway to determine the dissipativc effects in heavy ion collisions [33] -[36]. 

From the theoretical perspective, it should be possible to compute the 
shear and bulk viscosities directly from QCD at finite temperature. In prac- 
tice, this is extremely difficult because QCD is generally a strongly interacting 
theory. Calculations can and have been done at extremely high temperatures 
where perturbation theory, applied to quarks and gluons, can be used on ac- 
count of asymptotic freedom; see [37] for shear viscosity and [38] for bulk 
viscosity. At extremely low temperatures, perturbation theory can again be 
used because the matter consists only of a very dilute gas of pions, and low 
energy pion dynamics is well understood. See [39] for massive pions; for 
masslcss pions, see [39] for shear viscosity and [40] for bulk viscosity. There 
have been a variety of other kinetic theory calculations of the shear and/or 
bulk viscosities at low to moderate temperatures in the literature in recent 
years [41]-[45]; these usually include only elastic scattering of one or a few 
species of hadrons. 

In the intermediate region, which may be loosely defined as 100 < T < 
400 MeV, neither the low nor high temperature approach is accurate. A few 
lattice QCD simulations have used the Kubo formulae [46] to compute the 
shear [47, 48] and bulk [49] viscosities just above the critical temperature of 
pure gluon/glucball matter. However, accurate lattice QCD simulations of 
the properties of hadronic matter are extremely time consuming and the final 
results are still likely to be far in the future. The reason is that the lattice 
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spacing a must be small enough to describe the properties of an individual 
hadron (a < 0.05 fm) while the box size L must be large enough to contain 
many hadrons forming the dilute gas (L > 10 fm). Hence the number of 
spatial lattice sites should be at least 200 in each direction. An interesting 
alternative approach to the intermediate region is a model of classical, non- 
relativistic quasi-particles with color charges [50]. 

Our goal in this paper is to provide a theoretical framework in which to 
calculate the shear and bulk viscosities of hadronic matter. This framework 
has the following features. 

1. It is relativistic. 

2. It allows for an arbitrary number of hadron species. 

3. It allows for both elastic and inelastic collisions. 

4. It respects detailed balance. 

5. It allows for mean fields and temperature-dependent masses. 

6. The viscosities and the equation of state are mutually consistent in the 
sense that the same interactions are used to compute them all. 

Obviously some assumptions or approximations must be made for the theory 
to be applied in practice. The essential assumptions are that quasi-particles 
are well-defined and that the elementary interactions are local. Thus our 
proposed theoretical framework goes well beyond the classic works of [51] and 
[39] which, although relativistic, considered only elastic collisions in dilute 
gases. The inclusion of not only resonances, but especially inelastic collisions, 
mean fields, and temperature-dependent masses are essential for an accurate 
determination of the bulk viscosity [52]. 

The outline of this paper is as follows. In section 2 we recall the basics 
of the Boltzmann transport equation. In section 3 we derive the integral 
equations for the viscosities by using the Boltzmann equation. In section 
4 we show how the Landau-Lifshitz condition plays a crucial role for the 
bulk viscosity. In section 5 we work out formulas for the viscosities in the 
relaxation time approximation. In section 6 we generalize the previous results 
to include mean field or potential effects and their significance for the bulk 
viscosity. In section 7 we apply the framework to the linear a model with 
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massive pions. As expected, the ratio rj/s has a minimum and the ratio C,/s 
has a maximum near the rapid crossover transition, which is more pronounced 
for larger vacuum cr masses. We conclude in section 8. The reader not 
interested in mathematical details is referred to sections 7 and 8 and to the 
appendix where the main formulas are summarized. 

Since the asymmetry between matter and anti-matter in high energy nu- 
clear collisions at RHIC is very small, so are the baryon and electric charge 
chemical potentials. Therefore, thermal and electrical conductivity are ne- 
glected in this paper. There inclusion is straightforward but tedious, and 
work to include them is in progress. 

For case and clarity of presentation we will, for the most part, display 
formulas that include only 2 — )> 2 reactions, both elastic and inelastic, and 
formation 2 — > 1 and decay 1 — > 2 of resonances. At certain key points in the 
paper we simply write down formulas for the more general cases. In addition, 
in practice we generally use classical statistics. The lightest hadron for which 
this would be the most significant approximation is the pion, but even then 
the difference between Bose- Einstein and classical statistics have been shown 
to be inconsequential for the transport coefficients for zero chemical potential 
[53]. While some of the material here is not original, it is basic to developing 
the theoretical framework. We present it to make the paper self-contained 
and to define our notation. 



2 Boltzmann Equation 

The rate for 2 — > 2 processes, that is, the number of reactions per unit time 
per unit volume of the type a-\-h^c-\-d\s 

_ 1 /• d^Pa d^Pb d^Pc d^Pd I . ,| ,v|2 
l + 6aJ 2E„(27r)3 2i?,(27r)3 2E,(27r)3 2E,(2^)3l-^^"''l"'^^l 

x(27r)^5^ (Pa+p, -Pc- Pa) fah (l + (-l)'^=/c) (l + {-if'^h) (1) 

1 f d^Pa d^Pb d^pc d^pd 
= TT^y (2^ (2^0^ (2^0^ (2^0^^^"'^"^^ 

Xfafb (l + (-l)'*Vc) (l + (2) 
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whence 

Here is the spin of particle a, etc. and the factor 1/(1 + 5ab) takes into ac- 
count the possibihty that the incoming particles are identical. The amplitude 
M. is dimensionless. There are either Bose-enhancment or Pauli-suppression 
factors in the final state. In the rest frame of the system the single-particle 
distributions /„ are normalized such that 

j^fa{^,p,t) = na{:>^,t) (4) 

is the spatial density of particles of type a. In thermal equilibrium 

/a(x, P, t) = ^^Ea_^^)/T _ (_1)2.„ (^) 

where T is the temperature and jia is the chemical potential of the particle. 

It is obvious that the rate is a Lorentz scalar. 

The rate for decay processes, that is, the number of decays per unit time 
per unit volume of the type a — )■ c + d is 

r d^pa d^Pc d^pd . , . 
rate = / — — r^-— — — ^ , — -^r Allo c, a 



x(27r)V ip^-p^-p,) /„ (l + {-If^^u) (l + (-l)'^V.) 

d^Pa d^pc d^pd , 
(27r)3 (27r)3 (27r)3 



whence 

ItjaZrljc^rljd 
This A1 has dimension of energy. 

How do the W's relate to cross-sections and decay rates? The relation- 
ships for the cross-sections are as follows. 

dn* Q^ii^spl^J I 
<^cr 1 1 , . ^,9 /„s 
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Here s and t are the Mandelstam variables (one can tell from the context 
whether t represents a Mandelstam variable or time). The differential cross- 
section in the center-of- momentum frame is da/dfl* while da/dt is usually 
written as a function of the invariants s,t,u. Thus 

W{a,b\c,d) - ^ ' ^ ^ff ^' ;^(27r)^^^(p„+p,-pc-p,) . (9) 

This agrees with the literature on relativistic Boltzmann equations when 
only elastic collisions are considered since then Pgnai — Pfnitiai- alternative 
form, which may be more useful for cross-sections that are not isotropic in 
the center-of- momentum frame, is 

W{a,b\c,d) = P§ff^{27rfS'{pa+Pb-Pc-Pa) (10) 

where 

4s (p*„itiai)^ = (s - ml - ml)'^ - Amlml , 
4«(Pfinai)' = {s-ml-mlf-Amlml. (11) 

For the decay a ^ c + d consider the particle a at rest. It will decay 
according to the usual exponential law. 

—Jf - = -^a^c+d na{t) . (12) 

One computes that 

r„^e+d = ^|A^(a|c,rf)|2 (13) 

where 

4"^a (Pfinai)^ = i^l - ^1 - vfilf - Amlml (14) 

so that 

2 

W{a\c, d) = ^^L_r,^,+,(27r)^5^ (p„ -p,- pa) . (15) 

Now we consider the Boltzmann equation. Taking into account both gain 
and loss rates we write it as follows. 

dfa , f d^Pb d^Pc d^pa 



at ' " ^" (27r)3 (27r)3 (27r)3 
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+ 

cd 



^^iy(a, 6|c, (l + (-l)^^Ve) (l + (-l)'^V<i) } 

-iy(a|c, (l + i-ir^fc) (l + (-l)^^Vd) } 

- (l + } (16) 

Due to detailed balance on the microscopic level, energy conservation, 
and chemical equilibrium as represented by A*a + A*6 = A*c + for 2- body 
reactions and by /^a = Aic + f^d for 2-body decays, we find relations between 
the forward and backward going rates. 

{l + Sab)W{c,d\a,b) = {l + 6,d)W{a,b\c,d) (17) 
W{c,d\a) = {l + 6cd)W{a\c,d) (18) 



The Boltzmann equation then becomes 

dfa . 1 / d^Pb d^p^ d^Pdu/f h\ ^^ 

^ ^ ■ ^ S ^ ^ w (^ w 

x|/c/<i(l + (-l)'^Va) (l + i-ir^h) 
-fah{l + {-ir'fc) (l + (-l)'^''/.)} 

+E/^^^(«IM){/c/.(i + (-!)-/«) 



-fa (i + i-ir^fc) (l + (-l)^^V.) } 



-/aA (l + (-l)''-/c) } . (19) 

For classical statistics, where the Bose and Pauli factors are dropped, the 
equilibrium phase space distribution is 

/r(x,p,i)=e-(^-'^")/^. (20) 

The Boltzmann equation then shortens somewhat. 



+i:/||||^(ciM)a-/.A} (21) 

The generalization to arbitrary reactions {i} ^ {j} with n particles in 
the initial state and m particles in the final state is now clear. In obvious 
notation the Boltzmann equation is 



+ ^a-Vfa= y dP,W m\{j})F [/] , (22) 



dfa 



where the prime indicates that there is no integration over the momentum 
of a. There is a statistical factor for identical particles in the initial state 

S^Hml (23) 

i 

and products of Bose-Einstein and Fermi-Dirac distributions as appropriate 

^ [/] = n n {/. (1 + (-1)" fi) - (1 + (-1)^^' /.)} • (24) 

i j 
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3 Viscosities 



Now we use the Boltzmann equation to calculate the viscosities. We re- 
strict ourselves to zero chemical potentials. We assume that the system is 
in approximately local equilibrium, with local temperature T{x) and flow 
velocity U^{x). In the Landau-Lifshitz approach, U^{x) is the velocity of 
energy transport while in the Eckart approach, U^{x) would be the velocity 
of baryon number flow [54, 55]. However, the net baryon number, elec- 
tric charge, and all other conserved quantum numbers are taken to be zero. 
Therefore one cannot use the Eckart approach. Another consequence of all 
conserved quantum numbers being zero is that thermal conductivity has no 
meaning. 

The symmetric energy-momentum tensor is written as 

T^"' = -Pg^"" + wUW + AT'^^ (25) 

where P = P{T) is pressure, s = dP/dT is entropy density, e = —P + Ts 
is energy density, and w = Ts = P + e is enthalpy density. These are all 
measured in a frame in which the fluid is instantaneously at rest. The ATi^" 
is the dissipative part. It satisfies the condition 

U^AT"'' = (26) 

on account of the Landau-Lifshitz definition of fiow. The entropy current is 

= sU^ (27) 

and is conserved if dissipative terms are neglected. The most general form 
of AT'^'^ is given by 

AT'"' = r] (p'^U'' + D^'U^' + lA'^'dpUp) - (A^'^dpU^ . (28) 

Here 

A^" = U^U"" - g""" (29) 
is a projection tensor normal to C/'^, and 

D^ = d^- U^U^d^ (30) 
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is a derivative normal to f/^. The rj is the shear viscosity and the ( is the 
bulk viscosity. In the local rest frame of the fluid 



A"'^ = 

A^^' = S'^ (31) 

and 

Do = 

A = di. (32) 

In this frame 

d,s^ = ^ (d.W + d,U^ - |5^^V . U)' + |; (V • U)^ . (33) 

Non- decrease of entropy requires that both viscosities be non- negative. 

We are assuming that the interactions are well localized in space and 
time. We are assuming that they are, for practical purposes, point or con- 
tact interactions. Then the energy-momentum tensor is written as a sum of 
independent contributions. 

r-w=i:/|^f^/.(-.P) (34) 

Allow the system to be slightly out of equilibrium. This means that U^{x) is 
not constant in space and time, but that departures from local equilibrium 
are small. Then we can write 

Ux,p) = fl'^iU^pVT) [1 + 0a(x,p)] (35) 

and so 

^T'"" ^s^-^faHUc.pVT)Mx,p) (36) 

where \(f)a\ ^ 1. There is a constraint on (6a(a;,p) in that the Landau-Lifshitz 
condition (26) must be satisfied. It is customary and natural to use the 
same tensorial decomposition of (28) when expressing 0a(a;,p) as a function 
of space-time and momentum. 

0„ = -AadpUP + C;, (D'^U'' + D^Ui" + fA^^apt/^) (37) 
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Here Aa in general will depend on the scalar UaP""- The tensor C^^ could in 
principle be a linear combination of ^f^j/ and p^iPv However, the former gives 
zero contribution. Therefore we write C^j, = CaP^Pv where Ca will in general 
depend on the scalar UaP""- 

It is worthwhile emphasizing that the expansion of AT^*^ and 0^ in terms 
of the first order derivatives of the flow velocity is only an approximation. It 
is referred to as the first order dissipativc fiuid dynamics. Inclusion of second 
order derivatives goes under the names of Miiller and Israel and Stewart. 
The second order theory is under intense investigation due to its usefulness 
in describing high energy nuclear collisions where space-time gradients are 
not necessarily small. It is also worth emphasizing that these same quantities 
are zero for an equilibrated system in uniform fiow. 

It is now a straightforward matter to equate the two expressions for the 
dissipative part of the energy-momentum tensor. It is advantageous to work 
in the local rest frame of the fluid. 

C = 5i:/|f3^/r(^^./r)A„(Ej (38) 

How do we determine the and C^? The idea is to use the Boltz- 
mann equation where the term dfa/dt + • V/a is evaluated using the local 
equilibrium distribution fl^{UaP" /T). This is nonzero whenever the flow ve- 
locity is changing in space-time. It will act as a source for the collision term 
on the other side of the Bolztmann equation. With classical statistics the 
Boltzmann equation reads as follows. 



bed 



+ E / ^7C^/c'^W^Wa,&){0c-0a-06} (40) 

be (^^) i^^) 
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The first task is to compute the left-hand side of the Boltzmann equation. 
With /^^ = exp(-?7^p'^/T) we have 



(41) 



Using the conservation equations for energy and momentum, dJT^^ = 0, and 
entropy (since the viscous terms are neglected at this order), d^s^ = 0, we 
may deduce that 

(42) 



which has solution 



(43) 



where ,^ is a function of T which is undetermined by Eq. (42). It can 
be determined by substituting the above expression into the conservation 
equations. This gives ^ = —sjcy = —v^, where cy = de/dT = Tds/dT is 
the heat capacity per unit volume and — dPj de is the square of the sound 
velocity. Thus 



(44) 



After substituting in the structure of the 0's and grouping terms we get 
{dpU^) - ipi'U'' + D^U^" + lA^'^dpUp) = (45) 

where 



1 



3EaT 
d^Ph d^pc d^pd 



{p^U^f (l - 3^;,^) - mi 

W{a, 6|c, d) {A, + Aa-A,- A,} 



+ + (27r)3 (27r)3 (27r) 
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and 



be 

1 f d^ph d^Pc d^Pd .equ/^ U\ M 

X {Coi^cPc + CM - CaApl - aKK} • (47) 

The bulk viscosity is very small or zero in several limits. The first is the 
conformal limit, which means that the theory has no dimensional parameters, 
such as mass or intrinsic energy scale. Then P ~ and — 1/3, and so 
the first term on the right hand side of Eq. (46), the source term, vanishes 
and so do the The second is the nonrelativistic limit of a single species 
of particle. Then P ~ m^/^T^/^ exp(— m/T) and = T/m (plus corrections 
of higher order in T/m). Once again the first term on the right hand side of 
Eq. (46) vanishes (to lowest order in T/m) and so the bulk viscosity should 
be very small. These arguments do not apply to the shear viscosity since the 
source term does not involve the equation of state. 



4 Landau-Lifshitz Condition 

The equation (46) does not have a unique solution as it stands. For exam- 
ple, consider elastic scattering for just one type of particle. Starting with 

one solution A{E) we can generate an infinite number of other solutions by 
making the shift A{E) — )> A'{E) = A{E) — a — bE, where a and b are ar- 
bitrary constants. These constants are associated with particle conservation 
(a) and energy conservation (6). This has been noted in the literature before. 
It may be restated in more physical terms. To return a system to kinetic 
and chemical equilibrium after a change in volume, one might either change 
the number of particles while keeping the average energy per particle fixed, 
or one might change the average energy per particle while keeping the total 
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number of particles fixed. Now it is apparent that this ambiguity is associ- 
ated with the Landau Lifshitz condition (26), which is also sometimes called 
the condition of fit when solving (46) . 

Consider an arbitrary set of particle species and all possible reactions al- 
lowed by the symmetries. Make the shift Aa{Ea) — > A'^{Ea) = Aa{Ea) — da — 
bEa- The constant b must be the same for all species of particle. The con- 
stants Ua are just like chemical potentials; they satisfy the same relationships 
among themselves. Since we are restricting our considerations to systems 
with zero net quantum numbers, such as electric charge and baryon number, 
it is obvious that the all zero, just as all chemical potentials are zero. 

The constant b acts like an inverse temperature and is as yet undetermined. 

Suppose that we have a particular solution A^^'^ to (46); does it satisfy 
the Landau Lifshitz condition (26)? The general solution would be Aa{Ea) = 
A^^{Ea) — bEa- Using Eq. (36) the Landau Lifshitz condition for the A term 
is 

Here it is useful to know the contributions to the pressure, energy density, 
entropy density and heat capacity from a single species of particle. 



/ 



"^'^ Eafl\EjT) 



(27r)3 

r3 



CVa 



Now the coefficient b is determined in terms of integrals of the particular 
solutions. 

- ^ E / ^smEa/T)E^AriE^) . (50) 

If the particular solutions already happen to satisfy the Landau Lifshitz 
condition, then 6 = 0. Substitution of Aa{Ea) = /lP^'(£'a) — bEa into Eq. 
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(38), with b as determined above, gives an expression for the bulk viscosity. 



Notice that if the particular solutions happen to satisfy the Landau Lifshitz 
condition then Eq. (51) reduces to Eq. (38). 

There is no ambiguity with the Ca in the shear viscosity because of the 
tensorial structure of the integrand in Eq. (47). Physically the reason has 
to do with the fact that shear viscosity is associated with the response to 
changes in shape at fixed volume whereas bulk viscosity is associated with 
the response to changes in volume at fixed shape. 



5 Relaxation Time Approximation 

Consider the Boltzmann equation (21). Let us suppose that all species of 
particles for all values of momentum are in equilibrium except for species 
a with momentum p^. Replace all phase space distributions / with their 
equilibrium values f^"^ except for fa, which we allow to be out of equilibrium 
by a small amount. Thus we write = f^'^ + dfa- This is the momentum- 
dependent relaxation time approximation. We approximate the Boltzmann 
equation by 

^^"^""^^'^"^ + V„ • V/„(X, t, = -uMSfai^: Pa) (52) 

where 



+ (2^(2^^("l^''^) + 5^ (2^(2^^^^""'^)^^ ^'^^ 

is the frequency of interaction. The equilibration time is defined as 

ra{E)^u;\E). (54) 

The deviation 6 fa is related to the function 0^ defined in Eq. (35) by 

Sfa{x, p) = /^(a;, p)(p{x, p) . (55) 
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Therefore we can substitute Eqs. (37) and (44) into Eq. (52) to solve for the 
functions A^^^ and Ca, where C^'' = CaP'^Pa- 



3TEa 




(56) 



Ca{E,) 



2TEa 



(57) 



The viscosities are now readily calculated using these results. As usual, 
it is advantageous to work in the local rest frame of the fluid. 



C = ^ E/ [(l - 3-.') - ^'fi^iEa/T) . (58) 



In the relaxation time approximation one must calculate the momentum 
dependent relation time or rate from Eq. (53) and then substitute into 
the above expressions and perform a one dimensional integration. These are 
generalizations of the formulas given in [56] to an arbitrary number of species 
of particles with energy- dependent relaxation times. 

As a further approximation one may calculate a mean interaction fre- 
quency uJa and an associated mean relaxation time fa = i^a"^ via 



However, there is really no need to make this approximation unless one only 
desires a rough order of magnitude estimate. Depending on the dynamics, 
the relaxation time may be highly momentum dependent. In either case 
it is clear that the particles with the longest relaxation time dominate the 
viscosities, since these particles can transport energy and momentum over 
greater distances before interacting. 

Now comes a subtle point. What if we were to include weak interactions 
in our considerations? Clearly the relaxation times for weak interactions 




■ra{Ea)f:%EjT) 



(59) 




(60) 
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are orders of magnitude greater than the relaxation times for the strong 
interactions, so they would dominate the viscosities. The answer is that 
one must evaluate the actual physical conditions to which the viscous fluid 
equations are to be applied. For example, in high energy nuclear collisions 
the size of the system is on the order of 10 fm while the lifetime is of order 10 
fm/c. Any electrons, positrons or neutrinos that might be produced by this 
system will simply escape and not interact with any of the hadrons. They 
cannot transport energy and momentum to another part of the system. In 
addition, due to the weakness of the interaction very few of them will actually 
be produced. Thus the weak interactions are irrelevant in this situation. 
The environment in the early universe or supernovae will most likely require 
inclusion of the effects of the weak interactions since the length and time 
scales are so much greater. 



6 Mean Field or Potential Effects 

Most hadronic models of hot matter involve temperature dependent mean 
fields, such as the Walecka model or the linear sigma model. In the absence 
of chemical potentials there should be no condensation of vector or tensor 
fields, only scalar fields. Then the single particle energies have the form 
Ea — + fh1{T) where fha{T) is a temperature dependent effective mass 
that arises from the mean fields, in other words, potential energy effects. 
This affects the bulk viscosity, but not the shear viscosity, in several ways. 
The Boltzmann equation acquires an extra term on the left side of Eq. 

(16) 

+ |- • V, - W^Ea ■ j /„(x, t, p) = C[f] (61) 

where C[/] represents the collision, formation and decay terms. The extra 
term involves the force F = dp/dt = —VxE. To calculate the viscosities the 
left side is evaluated with the local equilibrium distribution 

f:^ = cxp[-U^{x)p:{x)/T{x)] . (62) 

Now p1 = Ea depends on x because m depends on T which depends on x. 
This approach is not new but has been proven or justified many times in the 
past. See, for example, references [57, 58]. 
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A straightforward calculation gives 

(- 



En 



1 



EnT 



J a 



eq 



dml 



dT 



The gradient of the temperature was determined previously. 



(63) 



(64) 



This allows us to write 



2EnT 



J a 



eq 



(65) 



Thus the term (1 — 3Vg)El — = — Sw^i^^ in Eq. (46) gets replaced by 
(1 - 3v',)El -fhl + ?,vlT^dfhl/dT^ = p2 - 3vl{El - T^dfhl/dT^). 

Taking into account the mean field effects modifies the particular solution 
in the relaxation time approximation. 



3TEa 



{l-3vl)El-ml + 3v',T' 



dT^ 



(66) 



However, the mean fields also affect the equation of state, the speed of sound, 
and the Landau-Lifshitz condition, so it is not so straightforward to deduce 
the bulk viscosity at this point. 



6.1 Model with no symmetry breaking 

Consider a model with N scalar fields that has no symmetry breaking, 
meaning that — for all a. It has the effective Lagrangian 

>Cefr = I E i^^^af - C/(^i, ^n) (67) 
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Assume that the potential is a polynomial in the fields to arbitrarily high 
order. It represents localized multi-particle interactions. The quasi-particle 
approach includes both mean fields and independent thermal fiuctuations 
around the mean fields, all calculated in a thermodynamically self-consistent 
manner [59]- [62]. This is sometimes referred to as the Phi- derivable ap- 
proach, and sometimes as the summation of daisy and super-daisy diagrams. 
If fermions are present they are integrated out and the effects of their inter- 
actions are subsumed in [/. In this approximation only even powers of the 
fields in the potential play a role. The thermal average may be written as 

(u)^ E ^m...n.($?r---«r (68) 

ni---njv 

where the Um-UN constants. The effective masses are obtained from 



a , 



where it is assumed that the system has been diagonalized in terms of normal 
modes such that 

^^S,,fnl (70) 



Prom the combinatorics 



/^\- V u (2n„)(2n.-l)(2n.-3)!! , 
\5$^/"„.^n/"™ (2n„-l)!!($^) ^^^^ ^^-^ ^''^ 

so that 



2 

m 



E ^n™n„($?)- • • • ($^)"- (72) 



\ al ni---njv 

The equation of state is given by 

P ^ Po-V 

e = eo + y (73) 

where the subscript refers to the free particle form with effective masses 
fha- The energy- momentum tensor is 

rj^f,^ = + gt^yy (74) 
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Note that Ts = e + P = eo + Pq = Tsq where s = dP/dT and Sq = dP^/dT 
so that the form of the entropy is unchanged. This is a consequence of the 
assumption of independent particle motion between coUisions in the mean 
field approximation. The potential energy density is obtained from 

a 

= E Un,-nM [1 - (ni + • • • + n^)] • • • {^ir- (75) 

ni---njv 

The entropy density is computed from the formulas 

dT V V dml dT dT ^ ' 



dPoa _ 

drhl ^ 



= -UK) (77) 



whereby it is readily shown that 



from which follows 



§ = E«o. (80) 



so that the model is thermodynamically consistent. There are some addi- 
tional interesting temperature derivatives that can be derived, such as 

(s.) 

The numbers Um-UN have the interpretation of representing the inter- 
action of 2{ni + ■ ■ ■ + un) particles in the initial plus final states. Taking 
account of the combinatorics, the vertex itself would be 

f^rti ■■■rtjv 

(2ni - 1)!! • • • {2nN - 1)!! 
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where it is assumed that all Ua 0. If an Ha happened to be zero then the 
factor (2na — 1) is replaced by 1. 

In Fermi liquid theory, functional variation of the energy density with 
respect to the distribution functions should yield the single particle energies 
[63]. That relationship holds here too. The kinetic part of the energy density 
is eo with the effective masses fha- Making a variation 6 fa also affects the 
mass. Hence 



d^p ME,) dml 5{^l) 
(27r)3 2E, d{^l) 5fa 



5V_ 



Variation of the potential part is given by 

5V = y. Un,...nA^-{ni + --- + nM)]{^lr 



(83) 



X 



ni 



(84) 



so that 

Sfa i^l) Sfa 



J2 Un,...n^ [1 - (m + • • • + tin)] • • • ($^)"- 

(85) 

Now it is easy to see from the expression (72) for the mass that dml/d{^l) 
is a symmetric matrix and that 



5V_ 



2 Sfa ^d{^l 

This cancels the extra term from eo so that 

— - p 



(86) 



(87) 



verifying the internal consistency of the model. 
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6.2 Model with symmetry breaking 

Now we allow for one of the fields to condense. For definiteness let it be the 
N'th one. After making the shift $jv '^'a^ + "^^j where v is the condensate, 
all fields obey ($a) = 0. Generalizing the previous analysis we write 

{U)= E Y.Un,...n.A^lr---{^lr''v' (88) 

ni---njv I 

If there is no condensation, then all Uni...nj^;i with / > vanish, and we return 
to the previous case. Following the same arguments as before we obtain the 
effective masses 

^a=7^ E Et^n,...n.;.n„($?r---($^)"V (89) 

and the potential energy density 

V = (C/)-iE^a(^a)' 

a 

= E E Un,-n.;i [1 - (ni + • • • + n^)] {^ir ■ ■ ■ {^ir^'v^ (90) 

The value of the condensate is determined by extremizing the pressure P — 
Po - ^ at fixed T, namely {dP/dv)T = 0. The result is 

E E Ur..-n.,il{^lr' ■ ■ ■ {^Ir-V' = (91) 
ni---nN I 

This is just the same condition as {dU/dv) = 0. 

For the dissipative part of the energy-momentum tensor and the Landau- 
Lifshitz condition we need 

= E E^n,...n.;ai-(ni + ---+n^)]($^)"^---«rV 
ni---njv / 

The condensate deviates from its equilibrium value because the thermal fluc- 
tuations deviate. One may express Av in terms of the A($^) using Eq. (91) 
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but for our purposes there is no need to do so explicitly. For verification of 
the Fermi liquid result that functional variation of the energy density with 
respect to the distribution function yields the single particle energy, it is 
sufficient to observe that 



This exactly cancels the extra term coming from eo so that 

JT = ^^ (94) 

in the presence of a condensate too. From this point on, the expression for 
AT^" and the bulk and shear viscosities is the same as when there is no 
condensation. 

6.3 Landau-Lifshitz condition 

The single particle energy at finite temperature is a functional of the dis- 
tribution functions Ea = Ea{{f}). When there is a small deviation from 
equilibrium 

fa {X,P)^ /r {E,,o) + Sfa {X, p) , (95) 

the single particle energy changes to 

Ea = Ea,o + SEa. (96) 

Here the subscript indicates the value the energy would have if there was 
no departure from equilibrium at all. If f^'^ is expressed as a function of the 
true energy Ea, then 

fa{x,p)^f:'^{Ea) + Sfa{x,p) , (97) 

where 

Sfa {X,p) = Sfa {X,p) - -I^^6Ea . (98) 

The structure of the quasi-particle Boltzmann equation is such that the it 
is the function 6 fa which determines the transport coefficients. It is impor- 
tant to realize that energy conservation in the collisions dictates that the 
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linearization of the collision term has to be performed with respect to the 
true energy Ea- 

Following the same arguments as earher, expansion of T*-' from Eq. (74) 
around local equilibrium, using Eq. (95), leads to 

The change in the energy density is given by 



AT' 



00 



Recall that 



and 



E 

o 

E 

a 

E 

a 

?/(27r)^ 
5fa = -e" 



{2nf 
(2^ 



EJfa 

(ejI + e^^I^^se., 



dEa 



E„ 



EaSfa - -ffr (Ea) SE, 



' ~ 1 dm^ ^ 

EaSfa-^^f!HEa)ST 



-Ea/T 



SEa Ea^ 



Sfa = e-^»/^ 



ST 

2^2 



Substituting these in (100) we find that, in the local rest frame. 



{2nfEa 

Hence, in a general frame of reference, 

d^p 1 



Ea - TTlaT-j^ ] Sfa 



dT 



{27rf Ea 



(100) 

(101) 
(102) 

(103) 



(104) 



This is the obvious generalization of the result of [58] to a system with mul- 
tiple species of particles. 



25 



Following the usual arguments we can deduce the viscosities. 



^4? 



V 



15 



(27r)3£;„ 



f:%EjT)A^iE,)\p\' 
fl\EjT)C,{E,) IpI^ 



(105) 



(106) 



(27r)3£;„- 

If one has a particular solution that does not happen to satisfy the 
Landau- Lifshitz condition, one can make it so by adding a term linear in 
the energy, just as before. 



E 



d^p 



(27r)3£;„ 



[Al^'{Ea)-bEa]=0. (107) 



To simplify the resulting formula, it is helpful to use 



E 



d^p 



f^iEjT) 



\p\^-3vl{El-T 



.dm^ 
dT2 



(108) 



which is a consequence of the identity dP/dT = {dP / de)de/ dT = v^de/dT in 
the independent particle models used here. The coefficient b is thus 



E 



d^p 



(27r)3£;„ 



fa^EjT) 



dT^ 



AriEa) 



(109) 



Substitution of Aa{Ea) = AP^{Ea) - bEa into Eq. (105) gives 



^4? 



d^p 



fl\EjT)Al^\E,) 



dT2 



(110) 



Of course, the term proportional to v1 will integrate to zero if the particular 
solution does satisfy the Landau-Lifshitz condition. The appearance of the 
factor in square brackets is natural since it matches the source function in 
Eq. (46). 

In the relaxation time approximation this becomes 



c = -E 



(27r)3 El 



f:\EjT) 



P - 3v? E^ - 



dT^ 



;iii) 



This expression is obviously positive definite. 
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7 The Linear a Model 



The linear a model has long been used as a simple renormalizable model of 
pion dynamics at low energy. Although it cannot claim to be quantitatively 
accurate, being supplanted by chiral perturbation theory, it is still a much- 
used model for testing approximations and as a proxy for more involved and 
detailed low energy models of QCD. Indeed, we will show that the model 
exemplifies the richness of the equation of state and transport coefficients to 
be expected from QCD. 
The Lagrangian is 

C=l{d,af + l{d,7rf-U{a,n) , (112) 

where 

U(a,n)^^{a' + 7r'-fy-Ha. (113) 

Compare to Eq. (67), where tTj is identified with $i for i — 1,2,3 and a 
is identified with $4. The SU{2)r x SU{2)l chiral symmetry is explicitly 
broken by the term Ha which gives the pion a mass. The scalar field has a 
non- vanishing vacuum expectation value v determined at the classical level 
by the equation 

\v(v''-f)=H. (114) 

The scalar field is thus split into a condensate and a fiuctuation, cr = v + A. 
The three parameters. A, H and /, are determined by the vacuum values of 
the pion decay constant ^-nd the pion and sigma masses. 

ml -ml 

A=^^ (115) 

H = flml (116) 

^ <fl^Sl (117) 

ml - ml 

For numerical calculations we take /tt = 93 MeV, m,r — 140 MeV, and either 
m^ = 600 or 900 MeV. 
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7.1 Thermodynamics 



The equation of state is calculated following the general procedure outlined 
in section 6. It must be done numerically. For more details the reader is 
referred to [59]-[62]. 

The temperature dependence of the meson masses and a condensate 
V are shown in Fig. 1 for two illustrative values of the vacuum a mass. 
There is no true phase transition, only a crossover from the low temperature 
regime where there is a large difference between the two masses and the 
high temperature regime where the masses are practically identical. The 
condensate decreases to very small values at high temperature but never 
vanishes. The temperature at which the symmetry is approximately restored 
is around 245 MeV. It should be noted for future reference that the pion 
and a masses rise linearly with temperature at high temperature. This is a 
generic feature of high temperature field theories. 

The entropy density, energy density and pressure are shown in Fig. 2. 
These thermodynamic quantities follow a continuous curve with no phase 
transition. They go to zero exponentially as T — )■ because all mesons are 
massive. At large temperature they show the usual behavior that s ~ 
and P ~ e ~ T'*. At intermediate temperatures rapid variation of the meson 
masses result in maxima in s/T^ and e/T^ but not in PjT^. 

Figure 3 shows the speed of sound squared v\ and heat capacity cy for the 
two choices of the vacuum a mass. As the vacuum a mass increases the speed 
of sound develops a dip and the heat capacity develops a peak around 245 
MeV. This indicates that the system is near a second order phase transition. 
At high temperature v\ — )■ 1/3 and cy — > T^, both on account of the fact 
that the equation of state approaches P T'^. 



7.2 Scattering amplitudes 

The Feynman rules can easily be determined from the Lagrangian of the 
linear a model. At one loop order the various scattering amplitudes are as 
follows: 

f Sabred + f^ac^bd H f^adhc , (US) 

s — t — u — J 
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Figure 1: (Color online) Variation of meson masses and condensate with tem- 
perature for a vacuum sigma mass of 600 MeV (top) and 900 MeV (bottom). 
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Figure 2: (Color online) Energy density, pressure and entropy density as 
functions of temperature for a vacuum sigma mass of 600 MeV (top) and 900 
MeV (bottom). 



30 




^ 4 



I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I 



I I I I I I 



TUa = 600 MeV ( 
rria = 900 MeV 




50 100 150 200 250 300 350 400 

T (MeV) 

Figure 3: (Color online) Temperature dependence of the speed of sound (top) 
and heat capacity (bottom). 
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M....a = -2A - + -^-^ + -^-^] , (119) 

\t — u — s — J 

M^^,,, = -2A - (^-^ + + -^-2] ' (120) 

\s — t — u — J 

Maa;aa = -6A - 36AV' f + + ) • (121) 

\s — t — u — J 

The poles in the s and u channels of tttt — > tttt and tto" — > ttct, respectively, 
pose trouble as is well known in the literature. We are not aware of any 
prescription to regularize both of these singularities while satisfying crossing 
symmetry. It is obvious that the problematic terms arc coming from the 3- 
point vertices. Those terms were not included when calculating the equation 
of state anyway. To be consistent with the equation of state, these terms must 
be dropped when calculating the viscosities. Equivalently we approximate 
the scattering amplitudes by their limits as s, t and u all go to infinity. Thus 
the scattering amplitudes just reduce to constants. 



7.3 Viscosities 

With the transition amplitiidcs specified, the departure functions Aa and Ca 
for a = TT and a can be determined as solutions to the integral equations 
Aa^O and C^" = using (46) and (47). When substituted into (105) and 
(106) we will get the bulk and shear viscosities. 

The set of integral equations are projected into a suitable vector space 
and solved by a variational method [51, 64]. The departure functions are 
expanded in terms of some chosen basis functions gn as 

N 

Ca^J2<'-,n9n{Ea/T) , (122) 
n=l 

and the expansion coefficients Ca^„ are adjusted to minimize the functional 

J'^J2{Ca-Ca,Ca-Ca) , (123) 

a 

where C is the true solution. We have defined the inner product between 
two functions Xi ^^'^ X2 as 

iXl, X2) = E / 7^Xl.a (P) ^ [X2,a (p)] ■ (124) 
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where 71 is the colhsion operator acting on the space of departure functions 
which, in the notation of (22), should be understood as 

1 „/ n m / \ 

nx]= E c/ dP^dP,w{{^}\{J})]J]Jfdl+{-lrfJ)[J:x^-T.x^ 

WO} i=i:'=i \ i j J 

(125) 

With the variational technique one hopes to find a good approximation to 
the departure functions using a modest set of basis functions. However, this 
is not guaranteed. In this paper we will show numerical results obtained 
this way for the shear viscosity. Solving the integral equations for the bulk 
viscosity is more difficult, as is well known, due to the presence of zero modes. 
Details of the numerical techniques used and their convergence properties will 
be presented elsewhere. 

In Figs. 4 {nia = 600 MeV) and 5 {rria = 900 MeV) we show the ratio 
of shear viscosity to entropy density from a 3rd order variational calculation 
compared to the relaxation time method. The variational method converges 
very quickly. For the present case the first order solution was already good to 
within a few percent and it was unnecessary to go beyond the 3rd order. It 
is satisfying that results from the relaxation time approximation are within 
a factor of 2 compared to the sophisticated variational method. A minimum 
near the crossover temperature of 245 MeV is seen in both cases as expected. 
The difference between the results of the variational method versus the re- 
laxation time approximation is most apparent at low temperature and for 
nicr — 900 MeV. This can be attributed to how the coUision dynamics is ap- 
proximated when using the relaxation time method. Although the relaxation 
time used is energy dependent, it is still only an approximation to the full 
solution of the integral equations. 
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Figure 4: Ratio of sliear viscosity to entropy density in tiie Cliapman-Enskog 
calculation (top) and in the relaxation time approximation (bottom). The 
vacuum sigma mass is 600 MeV. 
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Figure 5: Ratio of shear viscosity to entropy density in the Chapman-Enskog 
calculation (top) and in the relaxation time approximation (bottom). The 
vacuum sigma mass is 900 MeV. 
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Figure 6: (Color online) Ratio of bulk viscosity to entropy density in the 
relaxation time approximation. 

The ratio of bulk viscosity to entropy density from the relaxation time 
approximation 



c 



(27r)3 El 



f:\EjT) 



iPl 



2 / - 2 



1 2 



(126) 



is shown in Fig. 6. Although we have not yet managed to solve the integral 
equations for the bulk viscosity, we do not expect the difference between the 
two approaches to be significantly different from what we obtained for the 
shear viscosity. For a vacuum a mass of 900 MeV, there is a clear maximum 
in (/s near the crossover temperature. To ascertain its origin it is only 
necessary to examine Eq. (126). The ratio is quadratically proportional to 
the violation of conformality, either in terms of speed of sound 



3 s 
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or in terms of bulk transport mass 




adml _ d 



dT"^ d^' 



A comparison of Figs. 3 and 7 with Fig. 6 clearly shows that C/'S is largest 
when the violation of conformality is greatest. Note that for the system 
to be conformal it is only necessary that P ~ and that m ~ T; it is 
not necessary for the effective mass to vanish. Intuitively, the bulk viscosity 
is enhanced when it is easy to transfer energy between kinetic motion and 
internal degrees of freedom, such as resonances or heavier mass particle, 
mean fields, particle production or absorption, or effective masses that vary 
strongly with temperature. All of these are playing a role here to a greater 
or lesser extent. 

8 Conclusion 

In this paper we constructed a theoretical framework for the calculation of 
the viscosities of hot hadronic matter. The framework has the benefit that it 
is relativistic, it allows for an arbitrary number of hadron species, it allows 
for elastic and inelastic collisions and the formation and decay of resonances, 
it respects detailed balance, and it allows for temperature-dependent mean 
fields and temperature-dependent masses. It is a consistent theory in the 
sense that the same interactions that are used to calculate the equation of 
state are used to calculate the viscosities. The significance of this is that 
the bulk viscosity, in particular, depends very dramatically on the equation 
of state and on the quasi-particle masses. The essential assumption is that 
quasi-particles are good degrees of freedom. The basic resulting formulas 
have been assembled in the appendix for ease of reference. 

As a nontrivial application we applied this theoretical framework to the 
linear a model using the physical vacuum pion mass. As the vacuum a 
mass increases, the crossover transition become sharper, the speed of sound 
develops a dip, and the heat capacity develops a peak. This is well-known 
physics, but it has significance for the viscosities. As expected, the ratio r]/s 
has a minimum at the crossover temperature while the ratio C,/s develops a 
maximum there for large enough vacuum a masses. 
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Figure 7: (Color online) The quantity fh\dfh\l dT"^ which enters the calcu- 
lation of the bulk viscosity for a vacuum sigma mass of 600 MeV (top) and 
900 MeV (bottom). 
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As a byproduct of our theoretical framework we generalized the relax- 
ation time formulas for the shear and bulk viscosities of Gavin [56] to include 
an arbitrary number of species of hadrons with energy-dependent relation 
times and temperature-dependent effective masses. These formulas alone 
can be applied very usefully without much numerical work. 

The theoretical framework presented here can be applied to increas- 
ingly sophisticated models of hadronic interactions. Work on doing so is 
in progress. Also in progress is the extension to include chemical potentials. 
The latter may not be important at RHIC and LHC but they will be for 
heavy ion coUisions at FAIR. 
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Appendix 

In this appendix we gather some of the basic results of the paper. For more 
details see the text. 

The Boltzmann equation is 

^ + 1- • V. - v.K • V,] /a(x,i,p) = E ;^ fdP,dPjW({i}\{j})F[f] , 

(127) 

where the prime indicates that there is no integration over the momentum 
of a. There is a statistical factor for identical particles in the initial state 

^=n^*' (128) 

i 

and products of Bose-Einstein and Fermi-Dirac distributions as appropriate 

^ [/] = n n {fj (1 + (-1)^' f^) - (1 + (-i)^^' m ■ im 

i j 
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The rate function is 



wimm = \Mi{i}\{j})\\ 



(130) 



The departure of the phase space distributions are expressed in terms of (f>a 
as 

/„ = f:'^ (1 + , (131) 

which furthermore has the tensorial decompostion 

= -AadpU" + C;:, {^D^'U" + D^U" + lAi'^'dpUf) . (132) 
In terms of these the viscosities are 



(133) 



(134) 



The departure functions satisfy integral equations. For example, for two 
body reactions 
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(i^Pfe (i^Pc cP'Pd 



. + 6aJ (27r)3 (27r)3 (27r)3 
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1 f 

2E.T 1 + 7 



/,^'iiy(a,6|c,d){A„ + A- Ae- , (135) 



d^Pb d^Pc d^pd , I . 



bed 



+ 5aJ (27r)3(27r)3(27r)3- 



X {^pX + ap^Pfe - CcP^cP: - CrfP^^Pd} • (136) 

For resonance formation and decay and multiparticle reactions see the text. 
The dissipative part of the energy-momentum tensor can be written as 



(27r)'^„ 



dT2 



(137) 
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where 6 fa is defined in the text. In the relaxation time approximation the 
viscosities can be written as 



V 



15T 



1 




d^p |p|^ 



■ra{E^)mEjT) 



(138) 



and 



(27r)3 E^^ 




I 

(139) 



The effective masses are calculated self-consistently using the same potential 
U and same approximations as when calculating the equation of state, namely 
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